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1 Introduction 


Digital harmonic oscillators are standard devices used in sound synthesis, 
control systems, test and measurement equipment, heterodyning, frequency 
conversion, and many other DSP applications. A harmonic or sine wave 
oscillator generates a pure tone without any overtones. In additive synthesis, 
a number sine waves are added to generate timbre the Hammond organ is a 
prominent example of this principle. 


Methods for generating sine waves in the digital domain fall roughly into two 
categories: 


1. Recursive oscillators, where the actual output value is computed from 
previous states 


2. Direct evaluation as a function of a phase accumulator (= non band- 
limited ramp) 


The scheme presented in this note belongs to the first category. 


There are many algorithms for recursive oscillators, with different charac- 
teristics in terms of accuracy (especially at low frequencies), stability, com- 
putational complexity, etc.[1, 2]. In this note I present a new quadrature 
oscillator as the only hyperstable, equi-amplitude quadrature oscillator that 
I am aware of. 


Recursive oscillators perform best at constant frequency. Modulation re- 
quires coefficient(s) updating and may often result in amplitude change. For 


an equi-amplitude quadrature oscillator, the situation is less problematic. It 
is somewhat difficult to add phase modulation (it can be done for a quadra- 
ture oscillator, but it is not very efficient). 


2 Recursion Scheme 


Let u, and v, denote the oscillator outputs. The goal is to provide recursion 
equations for a harmonic quadrature oscillation, i.e. 


Un = cos(nw), Uy, = sin(nw) (1) 


with a specified frequency w in radians per sample. Given the state at some 
discrete point in time n, the next state may be found by a simple rotation. 
For reasons to become clear later, we will first advance from time n by half 
a sample to the intermediate time n + 1/2,: 


Un+1/2 = CUn — 8Un (2) 


Un pi/s = Sy CU, (3) 


where we have used the abreviations c = cos(w/2) and s = sin(w/2). Advance 
another half step to time n + 1, 


Un+1 = CUn+1/2 — SUn+1/2 (4) 


Un+1 = 8Un41/2 + CUn41/2 (5) 


Likewise, we may write down equations for the reverse transitions by half a 
sample frm n+ 1 ton + 1/2, 


Un+1/2 = CUnt1 + 8Un41 (6) 

Un+1/2 = —SUnti + CUn+1 (7) 
and from n + 1/2 to n: 

Un = CUn+1/2 + SUn+41/2 (8) 

Un = —SUn41/2 1 CUn41/2 (9) 


We will not make use of all eight equations (2)-(9), since they are redundant. 
It will be sufficient to use a selection of four, namely equations (2) and (6), 
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and equations (5) and (9). Thus, equation (2) provides an expression for the 
intermediate state tn41/2. Then equation (6) advances u another half step 
to the time n+ 1. To propagate v from n to n + 1, subtract equation (9) 
from equation (5) to obtain 


Un41 = Un + 28Un+1/2 (10) 


Note that the intermediate state un41/2 is only an internal temporary vari- 
able, a stepstone to leapfrog from v, to Uni. We may save two multiplies 
if we work with the auxiliary quantity w, := Un+i/2/c. Then equations (2), 
(10) and (6) become 


Wn = Un — kn (11.1) 
Unt1 = Un + koWn (11.2) 
Un+1 = Wn — kyUn41 C3) 


where the two constants k; and ke are 
ky =tan(w/2), ky = 2sin(w/2) cos(w/2) = sinw = 2k,/(1+k?) (12) 


Equations (11.1)-(11.3) constitute the main result of this paper. With initial 
values chosen as 


Ug = 1, Vo = 0 (13) 
the result of the recursion is a quadrature harmonic oscillation as in equation 
(1). 

Figure 1 shows the block diagram corresponding to equations (11). 
Below is a sample-by-sample processing algorithm in pseudocode. 


// initialize u and v 
at start dof 


u 1; 
v = 0; 
} 


// update coefficients 
if frequency changes do{ 
update w; 


ki = tan(0.5*w); 
k2 = 2*k1/(1 + ki*k1); 
} 


// iterate filter 
for every sample dof 
tmp = u - kikv; 

=v + k2*tmp; 
tmp - kixv; 


wes 
Il 


3 Stability 


For a stability analysis we take the z-transform of equations (11.1)-(11.3), 


w(z) = u(z) — kyu(z) 
zu(z) = v(z) + kow(z) 
zu(z) = w(z) — zkyv(z) (14) 


Eliminating w(z) and rearranging, we obtain the following homogeneous sys- 
tem of equations for u(z) and v(z): 


kou + [(1 = ky kg) = zu = 0 
(z —l)ut (24+ lk = 0. (15) 


Nontrivial solutions exist where the determinant vanishes, 


For 0 < kyky < 2, the roots 21,2 of equation (16) are complex conjugate and 
lie exactly on the unit circle. We may write 


a3= exw (17) 
with the frequency w given by 


cosw = 1 — kyko. (18) 


Even if k; and kg are subject to quantization in a real implementation, the 
poles still lie exactly on the unit circle hence preventing exponential runaway: 
the presented oscillator is hyperstable, as opposed to the well-known coupled 
form quadrature oscillator. 


Futrhermore, by virtue of equation (15), 


u/v = tisinw/ky = tiky/tan(w/2) = +i4/ ky (= — hs), (19) 


i.e. the eigenmodes exhibit a phase lag between u and v of exactly 7/2, 
regardless of coefficient quantization. The ratio of the corresponding ampli- 
tudes is unity if 

ky = 2k, /(1 + 7) (20) 


In practice, equation (20) will be valid only within machine precision, result- 
ing in slightly different amplitudes. 


4 Modulation 


An interesting question is how the system responds to frequency modulation, 
i.e. coefficient changes while running freely. Neglect for a moment roundoff 
errors introduced by finite precision arithmetics, then the trajectory in u-v- 
space will be an ellipse. A deviation from a circular trajectory occurs because 
equation (20) is fulfilled only to machine precision. Now if we change the 
values of k; and ko, the system’s trajectory will simply switch to another 
ellipse, possibly moving at a different speed. As long as equation (20) remains 
approximately valid, the new ellipse will be close to circular. 


If coefficients are changed frequently, the deviations may accumulate in the 
long run, ultimately resulting in an altered amplitude. Whether or not this 
is a concern depends on the application. For sound generation in a musical 
context, this is not an issue at all. 


5 Roundoff Error 


A rigorous analysis of the effect of roundoff errors on the iteration equation 
(11) is beyond the scope of this note. Quantization is often modelled as 
additive noise, without taking into account a possible correlation to the sig- 
nal itself. Such a model would result in a slow drift with the accumulated 
error increasing as the square root of the number of iterations. However, 
this behavior is not supported by empirical findings, which suggest that the 
system rather locks in into a limiting cycle. The oscillator has been run at 
44.1 kHz sample rate for days, with varying initial conditiona and coefficient 
values. In no case has there been an unbounded drift. Again, this is in harsh 
contrast with the coupled form oscillator, which shows systematic deviations 
already in the first few seconds, and ultimately terminates in an exponential 
runaway. 


Appendix: Some Common Oscillators 


Table 1 lists some recursive oscillators with their respective properties. For 
each listed oscillator tyxpe we give the governing equations below. 


Biquad Oscillator 


The biquad oscillator is a direct form I realization, and is the most used 
oscillator form. With one multiply and one add it is a very economic scheme, 
however accuracy is bad at low frequencies. 


Recurrence Parameter Initialization Output 
the = kt — tig RH 2cosw 1 = 0 Un = sin(nw) 
Up = — sinw 


Reinsch Oscillator 


This is another CPU friendly oscillator, with good accuracy at low frequen- 
cies. It is the ideal static oscillator. 


Recurrence Parameter Initialization Output 
Un+1 = Un + Un eS 4sin?(4w) uo = 0 Un = sin(nw) 
Oped Op — Pig gy Up = sinw Un, = Acos|(n + 5) 


A = 2sin(4w) 


Digital Waveguide Oscillator 


Another one-multiply oscillator with quadrature output, however with bad 
accuracy at low frequencies. 


Recurrence Parameter Initialization Output 


Sn = k(Un tn) k= cosw Up = 0 Un = Asin(nw) 
tigi = Ba Ue Uo = 1 Un = cos(nw) 
Un+1 = 8n — Un A=-— tan(5w) 

Unt1 = tn 


Quadrature Oscillator with Staggered Update 


Recurrence Parameter Initialization Output 

Unt] = Un + kup k = cosw Up = 0 Un = Asin(nw) 

Unti = kUn4i — Un vo = 1 Un = cos(nw) 
A=-—sinw 


Magic Circle Oscillator 


A popular ’almost quadrature’ oscillator. The outputs have equal amplitudes 
and a phase difference of 90 degrees plus half a sample. 


Recurrence Parameter Initialization Output 
Unt = Un — Kvn k = 2sin($w) Up = cos($w) Un = cos[(n — $)w] 
Unt1 = Un + kun vp = 0 Un = sin(nw) 


Coupled Form Oscillator 


The standard quadrature, equal amplitudes oscillator. Unfortunaltely it is 
not numerically stable: if left run freely, amplitudes will grow or decay ex- 
ponentially because the poles generally do not lie exactly on the unit circle. 
Therefore, some sort of automatic gain control has to be applied. Updating 
the parameters, i.e. frequency modulation, is somewhat CPU-demanding. 


Recurrence Parameters Initialization Output 


Unt = kiUn — kev, =3ki =cosw uw =1 Uo = cos(nw) 


Unt1 = kottn tkyvn kg =sinw wu=0 Up = sin(nw) 


Present Work 


A stable quadrature oscillator with equal amplitudes, good accuracy at low 
frequencies, and reasonable CPU load. 


Recurrence Parameters Initialization Output 
Wn = Un — k1Un ky = tan(5w) ug =1 Ug = cos(nw) 
Un+1 = Un + koWn ko = sinw v= 0 Uo = sin(nw) 


Un+1 = Wn — kyUn41 
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Figure 1: Block diagram of the oscillator corresponding to equations (11) 


Oscillator type Equal Quadra- Stable Low- 
ampli- ture frequency 
tudes accurate 

Biquad yes no yes no 

Reinisch no nearly yes yes 

Digital Waveguide no yes yes no 

Quad Staggered Update no yes yes no 

Magic Circle yes nearly yes yes 

Coupled Quad yes yes no yes 

This work yes yes yes yes 


Table 1: Some common oscillators and their properties 
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